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Abstract. We develop a mean-field model describing the Hamiltonian interaction of ultracold atoms and 
the optical field in a cavity. The Bose-Einstein condensate is properly defined by means of a grand- 
canonical approach. The model is efficient because only the relevant excitation modes are taken into 
account. However, the model goes beyond the two-mode subspace necessary to describe the self-organization 
quantum phase transition observed recently. We calculate all the second-order correlations of the coupled 
atom field and radiation field hybrid bosonic system, including the entanglement between the two types 
of fields. 
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1 Introduction 

A thermal cloud of cold atoms interacting with a single 
mode of a high-finesse optical cavity undergoes a phase 
transition when tuning the power of a laser field which 
illuminates the atoms from a direction perpendicular to 
the cavity axis. Below a threshold pump power, the cloud 
is homogeneous which is stabilized by thermal fluctua- 
tions. In this phase, the optical mean field in the cav- 
ity is zero, because the laser pump field is not scattered 
into the cavity mode from the homogeneous distribution 
of atoms. Above a threshold pump power, however, this 
solution becomes unstable. Then the atoms self-organize 
into a wavelength-periodic crystalline order which gives 
rise to Bragg-scattering from the transverse pump laser 
into the cavity. The resulting cavity field traps the atoms 
in the optical lattice distribution (see Figure 1). The self- 
organization effect has been first predicted in Ref . [I] , and 
soon experimentally observed [2]. It is closely related to 
the collective atomic recoil lasing transition [3] which is 
a Kuramoto- model- like synchronization phenomenon [4]. 
A classical mean-field description of self-organization has 
been presented in [5,6], which relics on a self-consistent 
canonical distribution of the atomic ensemble at a finite 
temperature T. The mean-field theory in [7] introduces a 
Vlasov-type equation in phase space which can account 
for arbitrary velocity distributions. 

The same self-organization effect can occur in the case 
of Bose-Einstein condensed ultra-cold atoms (BEC) at 
zero temperature [8,9,10,11,12]. For low pump power, in- 
stead of the thermal fluctuations, the homogeneous phase 
is stabilized by atom-atom collisions, or, in the lack of 
collisions, ultimately by the zero-point kinetic energy. A 
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sharp transition threshold is thus expected at zero tem- 
perature, too. This has been experimentally evidenced re- 
cently [13]. Moreover, in Ref. [13] as well as in Ref. [14], an 
analogy of the self-organization to the celebrated Dicke- 
model phase transition [15,16,17] has been pointed out. 
The dynamically coupled BEC wavefunction and single- 
mode cavity field realizes the Dicke-model with tunable 
parameters in the kHz range, and this is the first system 
in which the critical point can be reached and investigated. 




Fig. 1. Self-organization of a Bose-Einstein condensate 
in a cavity. Below a certain threshold pump power (left), 
the ultracold atoms have a quasi-homogeneous distribu- 
tion, and the cavity field with mode function cos(fcic) is 
empty. Above threshold (right), the atoms self-organize 
into a A-periodic ordered lattice in which they scatter the 
pump light constructively into the resonator mode. There 
is another, --shifted lattice possible to be formed. 



The Dicke-type phase transition is considered usually 
in systems with a fixed number of atoms and where the 
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individual atomic degrees of freedom span only a 2-mode 
Hilbert space. These assumptions are necessary to intro- 
duce the spin representation of the two-mode boson field. 
The drawback of this approach is that the role of higher 
excited motional modes cannot be included in the descrip- 
tion. 

In this paper we resort to a different approach which 
allows for the generalization to a multimode treatment 
of the matter wave field. Instead of fixing the number of 
atoms in the atomic modes, we invoke the grand canonical 
ensemble description, on assuming that the constant N c is 
the mean number of atoms in the condensate. In contrast 
to Ref. [8] where the condensate mean field has been writ- 
ten in position space, here both the condensate and the 
quantum fluctuations will be treated in momentum space. 
This is the most economic approach in terms of compu- 
tational needs, since the calculation converges very fast 
to the exact result as one includes higher excited kinetic 
energy eigenstates. We will show that the position of the 
critical point is not affected by the higher modes. Further- 
more we will confirm that the two mode model correctly 
describes the system and the phase transition below, and 
in the vicinity of the critical point. Far above threshold, 
however, the effect of the higher modes will become sig- 
nificant and the new multimode treatment is required. 

The paper is organized as follows. In section 2, we de- 
scribe the microscopic Hamiltonian model of our system. 
In section 3, the one-particle wavefunctions and the mode 
expansion are introduced. The backbone of our paper is 
section 4, where we introduce the grand canonical Hamil- 
tonian, which makes possible to systematically define, in 
subsection 4.1, the mean-field approximation in a multi- 
mode model. Then, in 4.2, we determine the independent 
quasiparticles by means of a Bogoliubov transformation. 
The subsection 4.3 is devoted to studying the fluctuations 
in the normal phase of the system (below threshold). In 
section 5, the ground state of the fluctuation Hamiltonian 
is analysed and the incoherent populations in the excited 
modes above the condensate are calculated numerically. 
We show that the ground state is an entangled one of the 
bipartite system of the cavity mode and the atomic mo- 
tional degrees of freedom, and the entanglement is quan- 
tified in section 6. 



2 Microscopic model of the system 

We consider an ensemble of ultracold atoms at T = in- 
teracting with a single-mode of a high-Q optical cavity 
[18,19]. The atoms are coherently driven from the side by 
a laser field with frequency u, directed perpendicularly to 
the cavity axis (see Fig. 1). The driving strength is de- 
scribed by the Rabi-frequency fin- The laser is detuned 
far below the atomic transition lja, that is, \Aa\ 2> 7, 
where 27 is the full atomic linewidth at half maximum 
and the (red) atom-pump detuning is Aa = uj — oja < 0. 
This condition ensures that the electronic excitation is ex- 
tremely low in the atoms, hence the spontaneous photon 
emission is suppressed. At the same time, the laser field 
is nearly resonant with the cavity mode frequency a>c, i.e. 



\Ac\ ~ k, where k is the cavity mode linewidth and the 
cavity-pump detuning is Aq = uj — ujc ( au these param- 
eters are summarized in Fig. 2). The scattering of laser 




laser 
co, Q R 



Fig. 2. Summary of the system parameters used in the 
microscopic model. 

photons into the cavity is thus a quasi-resonant process. 
Moreover, it is significantly enhanced by the strong dipole 
coupling between the atoms and the mode due to the small 
volume of the cavity. This coupling strength is character- 
ized by the single-photon Rabi frequency g, which is in 
the range of k. 

For the sake of simplicity, we describe the dynamics 
in one dimension x along the cavity axis, in which direc- 
tion the mode function is cos(kx), and the cavity length is 
L. The atom field and the resonator mode are described 
by the pair of bosonic annihilation and creation field op- 
erators, &(x), &l(x) and a, a\ respectively. In the large 
detuning limit, i.e., Aa is the far largest frequency in the 
system, the excited state can be eliminated [20], and the 
many-particle Hamiltonian in the reference frame rotating 
at the pump frequency u> is (h = 1) : 

+ Uq a cos 2 (fcx) + T] t cos(kx) (af + a)W(:r) dx . 

(2.1) 

The first term gives the optical field energy in the res- 
onator. The second term is the kinetic energy of the atoms. 
The third term describes the dispersive interaction be- 
tween the atoms and the cavity with a coupling strength 
Uq = g 1 1 Aa- The underlying physical process is the ab- 
sorption and stimulated emission of a cavity photon. This 
scattering process, from the point of view of the atoms, 
means an effective potential of the shape cos 2 (/cx) with 
depth depending on the photon number a) a. According to 
the red atomic detuning {A a < 0) the coefficient Uq < 0, 
which implies trapping positions at the A/2 separated 
antinodes of the cavity mode function. On the other hand, 
from the viewpoint of the cavity field, this term preserves 
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the photon number but gives rise to a frequency shift de- 
pending on the spatial distribution of the atom field. This 
term is responsible for the optomechanical-type coupling 
[21,22,23] investigated in the experiments [24,25]. The last 
term describes the effect of the pump field and results 
from the scattering between the pump laser and the cav- 
ity field. The back action of this scattering on the pump 
laser is neglected. It thus amounts effectively to a cavity 
field driving, + a, with the constant transverse pump 
amplitude r\ t = fljigjA^ and depending again on the 
local matter wave field density. 

In Eq. (2.1) we consider only one type of atom-atom 
interaction, namely the one which is mediated by the cav- 
ity field. This interaction is long-range [26] and has strict 
periodicity due to the fixed momentum of the exchanged 
photons. In contrast, the effect of s-wave scattering can 
take place at arbitrary momentum values and causes a 
broadening of the atomic momentum distribution around 
the peaks at integer times the momentum of the cavity 
photon [27]. We disregard this broadening effect by as- 
suming that the cavity mediated long range interaction 
dominates over s-wave scattering. Such an assumption is 
physically sound for a wide range of experimental para- 
maters [21] and helps the distillation of the effect caused 
by the cavity photons. 

Here we are interested in the ground state proper- 
ties and the excitation spectrum of this Hamiltonian, and 
therefore disregard the effects of the photon leakage out 
of the cavity. 



3 Mode expansion 

For the decomposition of ^f{x), we can use the complete 
orthonormal set of mode functions, 




cos(nkx) 




sin(nkx) 



However, the sin(nfc:r) modes are not populated by the 
parity-conserving Hamiltonian (2.1) when the system 
starts from a homogeneous BEC (being the ground state 
for Uq = rjt = 0). Thus, ^(x) can be expanded in terms of 
the even modes as 



E 

n=l 



— cos(nkx) c n 



(3.1) 



where the c n operators satisfy bosonic commutation 
relations. We will use the compact notations c = 

(co , c\ , . . .) T , being a column vector, and = ^cj , c\ , . . 

being a row vector. Then the operator of the total atom 
number is: 



N 



& (x)W(x) dx 



71=0 



c 



(3.2) 



On substituting the expansion (3.1) into (2.1), the Hamil- 
tonian can be constructed in terms of quadratic forms as 



H = - A n a f a - 



V2 



LOR 



M<°> c 



T] t (a f + a) (c f M (1) , 
- U a f a (c f (m (2) +21 



(3.3) 



where lur = ^— is the recoil frequency, I is the unit matrix, 
and the matrices are all real and symmetric, 



M<°> 



/0 2 



V 

/o i 

1 

1 

v/2 



(3.4a) 
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(3.4b) 
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10 1 
V20 01 
1 00 
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V ■ ■■) 



(3.4c) 



Note that the matrix is diagonal for j = 0, tridiag- 

onal for j = 1, and pentadiagonal for j — 2. Since the ki- 
netic energy difference between adjacent modes increases 
with the mode index, there appears a natural cutoff ex- 
cluding the high energy modes from the dynamics and the 
system effectively has only a finite number of atom field 
modes. In this way, this mode decomposition will result in 
a significantly reduced numerical effort compared to the 
real-space mean-field approaches in Refs. [8,28]. 



4 Mean-field approximation 

We assume that the atoms form a condensate which is a 
macroscopic mean-field (ip(r,t)) rotating at a unique fre- 
quency corresponding to the chemical potential (i. There- 
fore, we transform into the picture given by the grand 
canonical Hamiltonian 



K = H - fiN 



(4.1) 



which defines a dynamics such that the condensate mean- 
field is static. The chemical potential can be determined 
self-consistently, by assuming a fixed density of the con- 
densate atoms N c /L, 
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Let us separate the mean values of the operators a and 
c from the quantum fluctuations. The mean-field part is 
denoted by (a) = ^/W c a and (c) = ^fW c ~f, where the 
7 T 7 = f normalization condition fixes the number of con- 
densed atoms to be N c . Later we will justify that a and 
the elements of 7 can be chosen real. The y/Nc multipli- 
ers are included in order to make a and 7 constant in the 
thermodynamic limit (N c — > 00, L — !• 00). Let us displace 
the operators, 



a — » \J N c a + a , 



(4.2a) 
(4.2b) 



which is a canonical transformation. After the displace- 
ment (a) — and (c) = 0, so the new operators cor- 
respond to the quantum-fluctuations. Note that the dis- 
placement breaks the U(l) symmetry of the microscopic 
model, which is the invariance of the Hamiltonian with 
respect to the transformation >P(x) —> <P(x)e~ lv> with ar- 
bitrary phase ip. The expectation values of the total pho- 
ton and atom number operators can be expressed in the 
displaced frame as 



(a f a) N c a 2 + {a) a) 



(c f c) 



(c f c) 



(4.3a) 
(4.3b) 



That is, on top of the number of iV c condensate atoms, 
there is an incoherent population of atoms outside the 
condensate which appears due to the atom-photon inter- 
action. 

The terms in K after the displacement should be 
grouped according to the powers of a and c, 

£(0) ,£(!)+ £(2) ,£(3) ,£(4) 



K 



(4.4) 



Let us introduce new parameters which have constant 
value in the thermodynamic limit {N c —¥ 00, L 00): 



Ac 



1 



N C U 



1 



(4.5a) 

(4.5b) 
(4.5c) 



u = - N c U 

y = ^2N~Vt 
The zeroth order term of K is a c-numbcr, 

K {0) = N c ( -5 C a 2 + uj r (~f T M (0) 7) 



+ y a (~f T M (1) 7) + u a 2 ( 7 T M (2) 7) - fij , (4.6) 

so it is irrelevant to the dynamics, but gives the mean-field 
energy of the system. The first and second order terms 
read 

K (1) = s/W c (a f + a) ^(7) a+^y~f T M (1) 7 ^ 

+ ^N~ c (c f (M(a) - fi I) 7 + 7 T (M(a) - fil) c) , 

(4.7) 



and 

K {2) = 0(7) a) a + c f (M(a) - fit) c 

+ - (a f + a) (c f M'(a) 7 + j T M'(a) c 

The effective cavity resonance frequency is 
^(7) = -8c +W7 T M (2) 7, 



(4.8) 



(4.9a) 



and the cross-coupling of the motional modes via the cav- 
ity mean field is expressed by the matrix 

M(a)=w fi M (0) +i/ttM (1) +Ma 2 (M (2) +2l) , (4.9b) 

which is a real, symmetric matrix valued polynomial of 
a. M'(a) denotes the derivative of this polynomial with 
respect to a. 

It can be seen that is proportional to Nc , 

which implies that the third and fourth order terms dis- 
appear in the thermodynamic limit, so we disregard them. 
The mean-fields a and 7 are determined by the condition 
that, in the grand canonical Hamiltonian, the terms linear 
in the fluctuations a and c collected in have to vanish. 
The fluctuations around the mean values are described by 
the bilinear Hamiltonian K^ 2 \ 



4.1 The mean field solution 



The condition K^- 1 ' = leads to the system of equations: 



Q{i)a + ^y7 T M (1) 7 = 
M(a) 7 = fi 7 . 



(4.10a) 
(4.10b) 



These equations define a quasi eigenvalue problem: 7 is the 
eigenvector of the matrix M(a), and the smallest eigen- 
value is the chemical potential. Because the matrix is sym- 
metric, fi will have a real value. But the value of a, and 
so the matrix itself depends on 7 through the first equa- 
tion, which renders the problem to be nonlinear. It can be 
solved by iteration which, as a main virtue of the present 
approach, is a stable and fast numerical method. 

Let us make an iteration step starting from the initial 
value a = 0. Then, M(a = 0) — lur M^ -*, which is a diag- 
onal matrix and its smallest eigenvalue is fi = 0. The cor- 
responding normalized eigenvector is 7 = (1 , 0, 0, ...) T 
for which J?(7) = -5 C + and 7 T M (1) 7 = 0. This 
yields a = 0, which is then a trivial solution describing 
the normal phase of our system: the resonator contains 
no photons and the whole condensate is in the homoge- 
neous mode. This solution always exists but it becomes 
unstable above a certain threshold y cr it = s/—Sc wr- This 
critical point can be seen in Fig. 3 which plots a, the 
mean amplitude of the cavity mode divided by y/N c , as a 
function of the transverse pump amplitude y. This curve 
is calculated for various n cu t ff cutoff mode indexes, from 
"cutoff = 2, corresponding exactly to the Dicke-model [14], 
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Fig. 3. The coherent field amplitude a as a function of 
the pumping strength y for various cutoff mode numbers 
("cutoff = 2 . . . 10). The critical point y cr j t = 10 does not 
depend on n cuto s- The parameters: u)r — 1, 8c = —100, 
u = -20. 




Vl Vcrit 

Fig. 4. The distribution 7? of condensate atoms in the 
modes as a function of the pump strength. The phase 
transition occurs in the subspace spanned by the modes 
Co and c\, but well above the threshold other modes also 
get involved in the dynamics. The parameters: = 1, 
6 C = -100, u = -20, n cutoff = 10. 



still far below the one needed in the real-space description 
[8]. 

The role of higher-order modes is illustrated also in 
Fig. 5. There is another criticality in the system, of dif- 
ferent nature, which occurs when the effective mode fre- 
quency J7, depending on the condensate distribution 7 
as shown in Eq. (4.9a), becomes negative. In this regime 
there is no stable solution for the coupled atom field and 
cavity mode system. The effective mode frequency can be 
tuned by varying |u|, which leads to a divergence in a. 
This critical point depends on the cutoff mode index be- 
low n cu toff- In the two mode case the divergence occurs 
when u — 8c, however, the exact result, approached well 
with the cutoff rt cu toff = 10, is at much smaller |u| because 
the atoms in higher order modes are allowed to localize 
much better at the antinodes of the cavity mode function 
and yield a larger resonance shift. 
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Fig. 5. The effect of the phase shift term on a. The nu- 
merical simulation confirms that the effective cavity fre- 
quency 1^(7) tends to zero at the divergence point. The 
parameters: lor = 1, 8c = —100. The value of y is fixed 
just above the critical point, y — 11, y cr i t = 10. 



to "cutoff = 10. The case " cu toff = 10 is close to being ex- 
act, since the higher excited modes have j n ~ with four 
digit precision for n > 10. The distribution of the conden- 
sate atoms in the modes is shown in Fig. 4. Below thresh- 
old, only the homogeneous mode n = is populated. At 
the critical point, the population in the mode cos kx be- 
gins to grow abruptly from zero with a finite slope, while 
the higher mode populations start slowly with vanishing 
derivative. Therefore, the two-mode approximation holds 
in the vicinity of the critical point. Well above threshold 
the other modes get populated. In Fig. 3, above thresh- 
old, the significant dependence of the slope on the number 
of excited modes taken into account underlines the impor- 
tance of the multimode approach as contrasted to the two- 
mode model. However, the mode number " cu toff = 10 is 



4.2 The analysis of the fluctuations 

The values of a and 7 are known from the numerical solu- 
tion of (4.10) which also provides for the eigenvectors 
associated with the eigenvalues Xj of M(a). The matrix 
is symmetric and real, so the eigenvalues and eigenvectors 
are real. The eigenvectors form a complete orthonormal 
basis: v^' • v^ J ' = 8jj, We can arrange the eigenvalues 
in increasing order: the smallest one is A = /1, and the 
corresponding eigenvector is v^ ' = 7. 

Let us first decouple the atomic modes interacting via 
the second term of (4.8). The eigenvectors put into the 
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columns of a matrix. 



O = v<°> 



A 2 ) 



(4.11) 



define the orthogonal transformation, O t O = OO t = I, 
which leads to the bosonic modes 

b = T -c. (4.12) 

Inversely, 

c = O b = v (0) • b + v (1) • h + v (2) ■ b 2 + . . . , (4.13) 

Since v(°) = 7 the bo mode describes the fluctuations par- 
allel to the condensate. Subsequently the b%, b 2 , ■•• modes 
describe orthogonal excitations. 

The grand canonical Hamiltonian simplifies to 

00 

3=0 

00 

+ 2 £ 9i (« f + «) (b\ + h) > (4-14) 

3=0 

where g = O t • M'(a) 7. The frequency of the 60 mode is 
zero, since Ao = n- This zero-mode is the Goldstone mode 
resulting from the U{\) symmetry breaking imposed by 
the choice of real mean field 7. In the two-dimensional 
phase space of the Goldstone mode, the quadrature f = 
^(6q + &o) is parallel with the condensate and is a constant 
of motion (commutes with the above K^>). The orthogo- 
nal quadrature \fW c 4> = *(&o~ bo)/ 2 corresponds to phase 



fluctuations of the condensate 
of motion 

dr 



29 



it 



dt 2 



It obeys the equation 
+ 0) , (4.15) 



(« 



(4.16) 



Let At be a small time interval on the timescale of the vari- 
ation of the condensate phase. 4>(At) can then be well ap- 
proximated by the second-order Taylor-expansion which 
includes the above two time derivatives. The growth of 
the phase fluctuations is characterized by (4> 2 (At)). Since 
(f)(0) is not correlated neither with a(0) nor with at(0), the 
only non-vanishing term up to second order is 



(^(At)) = (^m 



_9o_ 
4N C 



((a\0) + a{0)f) ■ At 2 

+ 0(At a ) , (4.17) 



where ((fit(O) + a(0)) ) will be given later by (4.18) and 
(5.2). However, this expectation value is close to 1 except 
for a small vicinity of the critical point where it diverges. 
In the thermodynamic limit N c 00 the phase undergoes 



then a free expansion with characteristic time scale about 
EvyvJ which is the far longest time scale. Note also that, 
below threshold go = so the phase fluctuations are not 
growing at all. 

From now on, we neglect these fluctuations and com- 
pletely eliminate the dynamics of the mode bo- This ap- 
proximation renders the condensate to be a classical back- 
ground field similar to the external laser pump field which 
was described by the real parameter r\ t in the model. 
This step is equivalent to projecting the atomic exci- 
tation space to the one orthogonal to the condensate, 
c — > Cj_ = c — 7(7 T c), as described in [30]. 

In the following, we perform a Bogoliubov-transforma- 

tion on , at , b , b ^ in order to define the independent 

quasiparticle modes , d^, which combine excitations 

of the atom field and the electromagnetic field. Let us 
introduce the quadrature amplitudes 



1 



Xo 



(a) + a) 



Po = i y j (a f - a) 



Pj 



(4.18a) 

(4.18b) 
(4.18c) 

(4.18d) 



where j G {1, 2, 3, ...}. Note that the (xo , Po) quadratures 
are related to a, not to bo- The quadratures obey the usual 
canonical commutation relations: 



[£k ,pi]=i Ski , 



(4.19) 



and all other commutators vanish. Using vector notations, 

x = (x a , X! , x 2 , ■■■) T and p = (p , pi , p 2 , ■■■) T , the 
grand canonical Hamiltonian can be expressed in terms of 
the quadratures as 



£(2) 



1 



1 



P P 



x S x + const , 



(4.20) 



2 r ^ 2 

where c-numbers were omitted and the S kernel matrix is: 



S = 



9x 



(Q 2 

gx (Ai - nY 

92 

V • 



!)2 



(A 2 - iuY 



7 



(4.21) 



with the off-diagonal elements gj = gj ■ \J f2(\j — fi). It 
follows that the Bogoliubov transformation amounts to 
the diagonalization of the real symmetric matrix S. Let 
U be the orthogonal matrix comprising the eigenvectors 



of S as its colums (U T 
transformation, 



U = U • U J =1). The canonical 



x = U-X. 



UP 



(4.22) 
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leads to 



oo 



, , 2 Y 2 
U 3 X 3 



const , 



(4.23) 



3=0 



where un are the real eigenvalues of the matrix S. This 
Hamiltonian describes independent harmonic oscillators 
associated with bosonic quasiparticles with Wj eigenfre- 
quencies. The annihilation and creation operators of the 
quasiparticles are 



(4.24) 
(4.25) 

(4.26) 




where j £ {0, 1, 2, ...}, then the Hamiltonian is 



£(2) _ ^ cjj d] dj + const 

3=0 



The spectrum, the set of the eigenvalues uij, is plotted in 
Fig. 6 as a function of the pump strength parameter y. 
The calculation relies on a number of modes n cuto g = 10 
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Fig. 6. The u)j eigenfrequencies of the system. For y = 
the well-known ujj — j 2 ■ ujr spectrum of an ideal gas in 
a box is rendered. On increasing y, the lowest eigenvalue 
tends to zero at the phase transition point. The parame- 
ters are the same as in Fig. 4. 



which reproduces the exact result [8]. 



4.3 Fluctuations in the normal phase 

We apply the general results to the normal phase of the 
system described by a = and 7 = (1,0,0, ...) T . 
So, there is no coherent mean optical field in the res- 
onator, and the condensate is homogeneous. Then (4.9a) 

and (4.9b) gives = -6 C and M(a = 0) = lj r M {0 \ 



respectively. This matrix is already diagonal: we can read 
out that /i = and \ n — n 2 ■ lur. The orthogonal trans- 
formation in (4.12) is the trivial one: O = I and b = c. 
The vector of the coupling constants between a and c is: 
g = O t • M'(a)7 = y (0, 1, 0, 0, ...) T . This is an im- 
portant result: in the normal phase a is coupled only to 
ci and the coupling constant is g\ = y. This means that 
the two-mode model [14] is exact in the normal phase, and 
the Hamiltonian (4.14) is simply 



00 

-6 C a) a + u R ^2 n2 4 c n 



+ \y(^+a) (c\+c 1 ) . (4.27) 



The S matrix is 



81 



S = 



V ■ J/crit 
y ■ 2/crit U)% 



4<4 



V 



9u% 



(4.28) 



with the modified coupling constant gi = y ■ 2/crit, where 
2/crit = V— 5c is the critical point of the phase transi- 
tion in the two mode model as shown in [14]. Diagonaliza- 
tion of the first block leads to the non-trivial eigenvalues 



(4.29) 




R 2 

J/crit 



The frequency cj_ goes to zero at y = y C rit, which is then 
the phase transition point. It follows that the phase tran- 
sition occurs at the point y = y CI - lt even in the multimode 
model. 



5 Analysis of the ground state fluctuations in 
arbitrary phase 

Our next goal is to express the ground state of the Hamil- 
tonian (4.8) in terms of the Fock space of the operators 
a and c which have clear physical meaning. The ground 
state of the Hamiltonian (4.26) is simply the vacuum state 
of the d operators. However, since the Bogoliubov trans- 
formation in Eq. (4.22) mixes creation and annihilation 
operators, the ground state contains photonic and mo- 
tional excitations, moreover, it will be an entangled state. 
As follows from the Hamiltonian in (4.27), below threshold 
the ground state is the two-mode squeezed vacuum. 

We can use the Wigner-function to fully describe the 
ground state [31]. We make use of the fact that the 
Wigner-function associated with the ground state of a 
bilinear Hamiltonian is always a multivariate Gaussian 
distribution which, when centered at the origin, is fully 
determined by its covariance matrix [32]. To obtain the 
Wigner-function, we have to calculate then the symmet- 
rically ordered covariance matrix. 
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We can start from the correlations between the (X, P) 
quadratures pertaining to the independent quasi-particles, 
which are 



(X k Xi) = - — Ski 

ZUJk 



(PkPl) = ^s kl 



2 



({X k P l } s )=0, 



(5.1a) 

(5.1b) 
(5.1c) 



where {...} s denotes symmetric ordering. 

Now we apply the (4.22) transformation to determine 
the covariance matrix of the (x, p) quadratures associated 
with the photonic and the atomic motional excitations: 



[XkXl) 



1 1 



3=0 
oo 



(PkPi) = -u)j 

3=0 

({XkPljs) = 



(5.2a) 

(5.2b) 
(5.2c) 



0.03 



0.02 - 



0.01 




Fig. 7. The number of the incoherent photons (a* a) (solid 
line) and the expectation value of the total number of 
atoms outside the condensate (c' c) (dashed line) near the 
phase transition point. The parameters: wr = 1, 6c = 
-100, u = -20, j/crit = 10, n cuto g = 3. 

As the application of the above results, let us calculate 
the number of incoherent photons in the ground state. 
From (4.18) and (5.2) it follows that 



3=0 



n 



n 



(5.3) 



This expression is numerically evaluated and plotted as 
a function of the pump strength parameter y in Fig. 7. 
With a similar calculation, just performing the replace- 
ments —¥ k and Q — > X k — (J,, we can get the number of 
atoms in each of the b k modes 



(blk) 



OO 



\ E K 

3=0 



(5.4) 



where k G {1,2,3,...}. By summing these terms, we can 
get the total number of atoms outside the condensate: 



oo 

k=l 

1 OO OO 

4ee^ 



k=i j=o 



Xk 



At 



(5.5) 



which is plotted in Fig. 7 as a function of the coupling 
parameter y. Note that the index associated with c k runs 
from 0, while the index associated with bk runs from 1. 
Finally, the number of atoms in the c k modes, which are 
associated with the spatial harmonic functions, can be ob- 
tained by use of the transformation rule (4.13), the inverse 
of the formulae (4.18) and the covariance matrix (5.2): 



O nk O nl (b{ bi) 




O n k O n l Uki Uu 



yj{\ k - fJ.){Xi - fi) 



(X k 


- (J>) 


(Xi 




/ (Xi 




V (A fc 


-A 



(5.6) 



Starting from this formula, and using the orthogonality 
of the matrix O, the total number of atoms outside the 
condensate in (5.5) can be verified. 



6 Entanglement in the ground state 

Let us partition the system to the cavity mode and an- 
other part including all the motional modes of the atom 
field. The ground state of the system |^ g ), which is the 
vacuum of the quasiparticlcs defined by the operators dk 
(k = 0,1,...), exhibits bipartite entanglement between 
the radiation and the matter wave fields. This is sim- 
ilar to the entanglement occurring in the Dicke-model 
[33,34,35]. The entanglement can be simply measured by 
the Neumann- entropy of the cavity subsystem C, 



Sc. 



-Tr c (pc • lnpc) 



(6.1) 



where the reduced density matrix of the cavity mode is 
obtained by tracing the total density matrix p = iV'gXV'gl 
over the atomic subsystem A, i.e., pc = Trj^ (p). Another 
entanglement measure, which can be used, is the linear 
entropy, 

S c ,un = l-Tr c (pi) ■ (6.2) 

Both entropies have zero value precisely if pc is pure. For 
mixed states, they are positive, however, Sc can have an 
arbitrarily large value, while the linear entropy is bounded, 
Sc,lin < 1- 
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0.6 



1.4 



V I Vcrit 



= 1, 



Fig. 8. The bipartite entanglement measures, the 
Neumann-entropy Sq and the linear entropy Sc u n as a 
function of the coupling parameter y. Parameters: 
6c = -100, u = -20, y crit = 10, n max = 3. 

Let us introduce the parameter 



(at + , 



(at 



a) 



V^oXpS). (6.3) 



which leads to a simple form of the entanglement measures 
in the ground state of the total system [36,37], 



c 2 V 2 



C, /in 



x-i lt /x-i 



i 

i — 



(6.4a) 
(6.4b) 



It can be seen that both functions are strictly increasing 
with increasing \- It follows from the definition Eq. (6.3) 
and from the Heisenberg uncertainty relation that \ G 
[+1, +oo), and the minimum x = 1, that is, a pure ground 
state corresponds to a minimal uncertainty state of the 
field mode. Here we exploited the relation (xrj) = (po) = 0, 
so that the quadrature variances are Axq = (.Tq) and 

Ap<3 = y/ (Po)- The entropies are plotted in Fig. (8), and 
they exhibit singularity at the critical point, similarly to 
the incoherent photon and atom numbers in Fig. 7. 



7 Conclusion 

In this paper we presented a ground-canonical mean field 
theory for the coupled system of a single cavity mode 
and many motional modes of an ultracold atom field. We 
identified all the necessary approximations to separate the 
mean field from the fluctuations. The theory has been ap- 
plied to describe the self-organization phase transition [8] . 
We showed that below, in the normal phase, and in the 
vicinity of the critical point the two-mode approximation 
for the atomic motion holds, therefore the system shows 
indeed an analogy with the Dicke model of superradiant 
quantum phase transition [14]. Well above threshold the 
multi-mode expansion is needed, however, our approach 



leads to mean-field equations which have to be solved nu- 
merically. We calculated the number of incoherent photons 
and populations in the higher excited motional modes, as 
well as the amount of entanglement between the matter 
and radiation fields. This approach is suitable and will 
be used in the future to deal with other type of multi- 
mode systems, for example, the matter wave field coupled 
to the radiation field in a degenerate confocal resonator 
which shows a rich phase diagram [38,39]. 
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